Open Tabs
- arbeidskrav 5
- kka120@ad.uit.no@jup01:~
- arbeidskrav 5
- Untitled.ipynb
- 5 - sympy.ipynb
- Untitled2.ipynb
Kernels
- Untitled2.ipynb
- Untitled.ipynb
- 5 - sympy.ipynb
- arbeidskrav 5
Terminals
- terminals/1
- assets2 months ago
- dataa month ago
- res2 months ago
- 1 - intro.ipynba month ago
- 1 - variabler.ipynb2 months ago
- 2 - funksjoner.ipynba month ago
- 3 - matplotlib.ipynb8 days ago
- 4 - lister, oppslag og numpy.ipynb8 days ago
- 5 - sympy.ipynb8 days ago
- 6 - logikk-filtrering-betingelser.ipynba month ago
- 6 - oppslag og pandas.ipynba month ago
- 7 - løkker.ipynb2 months ago
- 7 - simulering.ipynb2 months ago
- arbeidskrav 54 days ago
- arbeidskrav 5.ipynb5 minutes ago
- tips_and_tricks.md2 months ago
- Untitled.ipynba month ago
- Untitled2.ipynb8 days ago
import matplotlib.python as plt--------------------------------------------------------------------------- ModuleNotFoundError Traceback (most recent call last) /tmp/ipykernel_207225/3527871001.py in <module> ----> 1 import matplotlib.python as plt 2 import numpy as np 3 import sympy as sy 4 5 def S(q): ModuleNotFoundError: No module named 'matplotlib.python'
q = sy.symbol('q')xxxxxxxxxx5 - Sympy¶
xxxxxxxxxxSympy er en veldig nyttig pakke innenfor samfunnsøknomi. Med den kan vi regne analytisk, det vil si med symboler. I samfunnsøknomi bruker vi mye matematikk som denne pakken kan hjelpe oss med. Vi starter med et eksempel på utregning av profittmaksimum:
xxxxxxxxxxEksempel med optimal mengde arbeidskraft¶
xxxxxxxxxxAnta at du er bedriftsleder, og lurer på hvor mange medarbeidere du skal ansette. For en gitt mengde arbeidskraft produseres
Eksempel 1:¶
return 60*L**axxxxxxxxxxa definerer produktiviteten til de ansatte. Desto høyere a er, desto mer produktive er de ansatte. Vi kan plotte denne funksjonen:
Eksempel 2:¶
import matplotlib.pyplot as pltText(0, 0.5, 'Produksjon')
xxxxxxxxxxSom vi ser er produktiviteten avtakende, kurven stiger mindre utover i diagrammet. Det er fordi a er mindre enn én. (Forsøk med andre tall for a!)
Når produktfunksjonen er definert, kan vi definere fortjenesten til bedriften, eller "profittfunksjonen" som vi gjerne kaller det i samfunnsøkonomi. La oss si at bedriften betaler hver arbeider w tusen kroner, og at faste utgifter er K.
Fortenesten er pris ganger solgt mengde, p*f(L,a) og kostnadene er lønnskostnader w*L og faste utgifter K. Profittfunksjonen blir da
Eksempel 3:¶
return p*f(L,a)-w*L-KxxxxxxxxxxLa oss se på den grafisk. Dersom prisen per enhet er hundre kroner, lønna er 2 500 per dag og bedriften har 300 000 i faste utgifter per dag, så ser profittfunksjonen slik ut:
Eksempel 4:¶
plt.plot(x,profit(x,0.9,2500,100,300000),label='Fortjeneste')<matplotlib.legend.Legend at 0x1a422af56c8>
xxxxxxxxxxVi bruker metoden med fig,ax=plt.subplots(), siden vi skal bygge på denne grafen.
Vi skal nå begynne å bruke sympy. Det første vi må gjøre er å definere hvilke symboler som skal behandles analytisk (altså som symboler og ikke flyttall). Det gjør vi slik:
Eksempel 5:¶
L,a,w,p,K=sp.symbols("L a w p K")xxxxxxxxxxMed symbolene definert, vil nå vår profittfunksjon vises analytisk når vi bruker de definere symbolene:
Eksempel 6:¶
profit(L,a,w,p,K)xxxxxxxxxxVi ønsker å finne ut for hvilken arbeidskraft fortjenesten er størst. Det er det høyeste punktet i figuren fra Eksempel 3. På dette punktet har profittfunksjonen ingen stigning, slik at den deriverte er null. For å finne dette punktet må vi først finne den deriverte. Det gjør vi i sympy med funksjonen diff(). Den tar to argumenter. Det første er funksjonen, det andre er den variabel vi ønsker å derivere med hensyn til.
Som vi ser av figuren i Eksempel 3, så har vi arbeidskraft L langs x-aksen, så det er denne variabelen vi ønsker å derivere med hensyn til. Den deriverte av profitt() med hensyn til L er dermed:
Eksempel 7:¶
d_profitt=sp.diff(profit(L,a,w,p,K),L)xxxxxxxxxxLøser førsteordensbetingelsen¶
xxxxxxxxxxFor å finne punktet der denne deriverte er null, setter vi opp en ligning der den deriverte er null, og løser for den L som tilfredstiller ligningen. En slik ligning kalles "førsteordensbetingelse", eller "first order condition" på engelsk. Vi forkorter den derfor til "foc":
Eksempel 8:¶
foc=sp.Eq(d_profitt,0)xxxxxxxxxxVi kan nå løse førsteordensbetingelsen med funksjonen solve, som ligger i modulen solversi sympy:
Eksempel 9:¶
from sympy.solvers import solvexxxxxxxxxxLegg merke til at resultatet ligger i en liste med lengde 1, så vi må hente ut element 0 i listen for å vise resultatet. Vi kan finne hva den analytiske verdien er i maksimum ved å sette L_maxinn i profittfunksjonen:
Eksempel 10:¶
profit_max=profit(L_max,a,w,p,K)xxxxxxxxxxVi kan nå beregne de nummeriske verdiene ved å sette inn noen passende tall for de ukjente symbolene. Vi prøver med 0.5 for produktivitet a, 0.3 for lønn w, 1 for prise p og 1 for kapital K. Vi forteller sympy om at vi ønsker å bruke dise verdiene ved å lage et oppslag der hvert av symbolene er nøkkel til hver av de nummeriske verdiene:
Eksempel 11:¶
num_dict={a:0.9,w:2500,p:100,K:300000}xxxxxxxxxxVi kan nå finne nummerisk hvor mye arbeidskraft som trengs for å oppnå maksimal fortjeneste:
Eksempel 12:¶
L_max.subs(num_dict)xxxxxxxxxxOm vi nå legger dette tallet inn for symbol L, kan vi finne hvor stor fortjenesten er i dette punktet. Vi starter med å legge inn verdien for L som gir maksimal fortjeneste:
Eksempel 13:¶
num_dict[L]=L_max.subs(num_dict){a: 0.9, w: 2500, p: 100, K: 300000, L: 2210.73919720734}xxxxxxxxxxMed det oppdaterte oppslaget blir fortjenesten
Eksempel 14:¶
profit_max_num=profit(L,a,w,p,K).subs(num_dict)xxxxxxxxxxMed modulene displayog Markdownfra IPython (pakken som driver Jupyter), kan vi sette dette in i en pen tabell:
Eksempel 15:¶
| Optimal mengde arbeidskraft: | ${np.round(float(num_dict[L]),1)}$ |${sp.latex(L_max)}$ | | Desimalverdi | Analytisk verdi | |
|---|---|---|
| Optimal mengde arbeidskraft: | ||
| Maksimal profitt |
xxxxxxxxxxLa oss nå plotte løstningen. Vi plotter her følgende, i rekkefølge:
- Den opprinnelige profittfunksjonen
- Den horisontale tangenten som tangerer i maksimumspunktet
- En vertikal linje som viser mengden arbeidskraft i optimum.
Legg merke til at vi legger inn verdiene vi har definert i num_dict inn i profittfunksjonen. Generelt er det en god idé i programmering å "hardkode" tall minst mulig. Definer det heller tallet som en variabel eller element i en dictog referer til det senere.
Eksempel 16:¶
ax.legend(loc='lower right',frameon=False)xxxxxxxxxxEksempel med tilbud og etterspørsel¶
xxxxxxxxxxI forelesning 3 definerte vi disse tilbuds og etterspørselsfunksjonene:
Eksempel 17:¶
return (x**2)*(1/250)xxxxxxxxxxOg vi tegnet dem slik:
Eksempel 18:¶
ax.plot(q,demand(q),color='green',label='Etterspørsel')xxxxxxxxxxVi løste da likevekten grafisk, ved å se sånn cirka hvor tilbud er lik etterspørsel. Med sympy kan vi la python regne ut dette, og konsument og produsentoverskudd. Vi gjør dette ved å definere mengde x som en eksogen variabel, og sette opp ligningen vi trenger, altså at tilbud skal være lik etterspørsel:
Eksempel 19:¶
eq_cond=sp.Eq(demand(x),supply(x))xxxxxxxxxxVi kan nå løse dette med solve fra sympy, som i forrige eksempel:
Eksempel 20:¶
x_eq=solve(eq_cond,x)xxxxxxxxxxBare én av disse løsningene er gyldige. De to siste i listen x_eqer såkalte "imaginære tall", det ser vi av I'en. Vi går ikke lenger inn på hva dette er her, men nøyer oss med å si at en likevekt ikke kan være et imaginært tall. Løsningen er altså x_eq[0]. Vi kan sette denne inn i enten tilbuds eller etterspørselfunksjonen for å få likevektsprisen
Eksempel 21:¶
Likevektskvantum er {x_eq[0]}xxxxxxxxxxEtterspørselskurven kan ses på som en rekke med konsumenter med ulik betalingsvilje i fallende rekkefølge. Alle konsumentene som betaler p_eq har dermed et overskudd som er lik differansen mellom p_eq og konsumentens punkt på etterspørselskurven. Summen av alle konsumentenes overskudd kalles konsumentoverskuddet. Dette kan illustreres ved å legge et skravert område til figuren over
Eksempel 22:¶
ax.fill_between(q,float(p_eq),demand(q), color = "pink",alpha = 0.3,label='Konsumentoverskudd')xxxxxxxxxxAkkrat som at vi kan regnet ut skjæringspunktet med sympy, så kan vi regne ut det skraverte konsumentoverskuddet. Vi bruker da det vi har lært i matematikkurset; arealet under en funksjon er integralet til funksjonen. Vi skal finne arealet under etterspørselsfunksjonen demand(x), men kun ned til prisen p_eq, så vi integrer differansen demand(x)-p_eq. Dette gjør vi for alle omsatte enheter, altså frem til omsatt kvantum x_eq[0].
Vi skal altså integrere demand(x)-p_eq i intervalet 0 til x_eq[0]. Det kan vi gjøre i sympy slik:
Eksempel 23:¶
consumer_surplus=sp.integrate(demand(x)-float(p_eq),(x,0,x_eq[0]))xxxxxxxxxxPå samme måte er produsentoverskuddet arealet over tilbuskurven, opp til prisen, altså det gule området i figuren under
Eksempel 24:¶
ax.fill_between(q,supply(q),float(p_eq), color = "yellow",alpha = 0.3,label='Produsentoverskudd')xxxxxxxxxxVi kan regne ut dette område også, som altså er integralet av differansen mellom prisen og tilbudskruven, frem til x_eq[0].
Eksempel 25:¶
producer_surplus=sp.integrate(p_eq-supply(x),(x,0,x_eq[0]))xxxxxxxxxxSummen av produsentoverskuddet og konsumentoversdkuddet kalles "velferdsgevinsten". Vi kan finne den ved å legge sammen konsument- og produsentoverskudd:
Eksempel 26:¶
producer_surpluss+consumer_surplussxxxxxxxxxxEller ved å ta integralet av differansen mellom etterspørsel og tilbud:
Eksempel 27:¶
welfare_surplus=sp.integrate(demand(x)-supply(x),(x,0,x_eq[0]))xxxxxxxxxxVi kan nå lage en tabell som oppsumerer resultatene:
Eksempel 28:¶
| Konsumentoverskudd: | ${np.round(float(consumer_surpluss),1)}$ | xxxxxxxxxxOppgaver¶
xxxxxxxxxxLøs oppgavene 15, 16, 17 og 19-23 i "Læresteg" i BED-1007 (se canvas for BED-1007), ved bruk av Sympy
xxxxxxxxxxoverskrift¶
xxxxxxxxxxfor counter in range(10)
print("%d er partall"%(counter))2 er partall 4 er partall 6 er partall 8 er partall
print("Please insert the current value: ")xxxxxxxxxx#oppgave 15.1a) import sympy as spx,y=sp.symbols('x y')eq1=sp.Eq(2*x + 4*y, 6)eq1eq2=sp.Eq(-2*x + y, 4)eq2from sympy.solvers import solvesol=solve([eq1, eq2], [x, y])sol{x: -1, y: 2}eq1.subs(sol)eq2.subs(sol)#beq1=sp.Eq(2*x - y, -1)eq1eq2=sp.Eq(x**2 + x - y, 1)eq2sol=solve([eq1, eq2], [x, y])sol[(-1, -1), (2, 5)]
eq1.subs(sol)eq2.subs(sol)#2a)eq1=sp.Eq(3*x - 12, -6*y)eq1eq2=sp.Eq(4*x - 8*y, 16)eq2sol=solve([eq1, eq2], [x, y])sol{x: 4, y: 0}eq1.subs(sol)eq2.subs(sol)#beq1=sp.Eq(x**2 + 4*x - 3, 3*y)eq1eq2=sp.Eq(2*y, 8 - x)eq2sol=solve([eq1, eq2], [x, y])sol[(-15/2, 31/4), (2, 3)]
#oppgave 16 1a)eq1=sp.Eq(4*x + 2*y, 12)eq1eq2=sp.Eq(6*x - 2*y, 8)eq2sol=solve([eq1, eq2], [x, y])sol{x: 2, y: 2}eq1.subs(sol)eq2.subs(sol)#beq1=sp.Eq(2*x + y**2, 25)eq1eq2=sp.Eq(x - 2*y, 10)eq2sol=solve([eq1, eq2], [x, y])sol[(0, -5), (12, 1)]
#2aeq1=sp.Eq(4*x - 4*y, 8)eq1eq2=sp.Eq(9*x - 5*y, 26)eq2sol=solve([eq1, eq2], [x, y])sol{x: 4, y: 2}#beq1=sp.Eq(x + 2*y**2, 15)eq1eq2=sp.Eq(x**2 - 4*y**2, 33)eq2sol=solve([eq1, eq2], [x, y])sol[(-9, -2*sqrt(3)), (-9, 2*sqrt(3)), (7, -2), (7, 2)]
#oppgave 17 1a)eq1=sp.Eq(x*y**2 - x, 0)eq1eq2=sp.Eq(x + y**2, 16)eq2sol=solve([eq1, eq2], [x, y])sol[(0, -4), (0, 4), (15, -1), (15, 1)]
eq1.subs(sol)eq2.subs(sol)#beq1=sp.Eq(x**2 + y**2, 100)eq1eq2=sp.Eq(x**2*y - 36*y, 0)eq2sol=solve([eq1, eq2], [x, y])sol[(-10, 0), (-6, -8), (-6, 8), (6, -8), (6, 8), (10, 0)]
#2aeq1=sp.Eq(x*y**2 - 49*x, 0)eq1eq2=sp.Eq(x**2 + y**2, 50)eq2sol=solve([eq1, eq2], [x, y])sol[(0, -4), (0, 4), (15, -1), (15, 1)]
#beq1=sp.Eq(x**2 + y**2, 5/4)eq1eq2=sp.Eq(2*x*y + y, 0)eq2sol=solve([eq1, eq2], [x, y])sol[(-1.11803398874989, 0.0), (-0.500000000000000, -1.00000000000000), (-0.500000000000000, 1.00000000000000), (1.11803398874989, 0.0)]
eq1.subs(sol)eq2.subs(sol)- arbeidskrav 5
- kka120@ad.uit.no@jup01:~
- arbeidskrav 5
- Untitled.ipynb
- 5 - sympy.ipynb
- Untitled2.ipynb
xxxxxxxxxx{ "cells": [ { "cell_type": "code", "execution_count": 11, "id": "3604f052-84c8-4cb0-b80f-e26495cb75df", "metadata": {}, "outputs": [], "source": [ "#oppgave 15.1a) \n", "import sympy as sp\n", "x,y=sp.symbols('x y')" ] }, { "cell_type": "code", "execution_count": 12, "id": "d4f92c29-c3a7-4856-aca4-6eca8e4ef245", "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$\\displaystyle 2 x + 4 y = 6$" ], "text/plain": [ "Eq(2*x + 4*y, 6)" ] }, "execution_count": 12, "metadata": {}, "output_type": "execute_result" } ], "source": [ "eq1=sp.Eq(2*x + 4*y, 6)\n", "eq1" ] }, { "cell_type": "code", "execution_count": 16, "id": "e37e2aa5-9aff-41c6-b6dd-b2826caea08a", "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$\\displaystyle - 2 x + y = 4$" ], "text/plain": [xxxxxxxxxxxxxxxxxxxxxxxxxxxxxx-
Variables
Callstack
Breakpoints
Source
xxxxxxxxxx